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Abstract. The Rayleigh-Taylor instability plays an important role in the dynamics of several astronomical objects, in particular, 
in supernovae (SN) evolution. In this paper we develop an analytical approach to study the stability analysis of spherical 
expansion of the SN ejecta by using a special transformation in the co-moving coordinate frame. We first study a non-stationary 
spherical expansion of a gas shell under the pressure of a central source. Then we analyze its stability with respect to a no radial, 
non spherically symmetric perturbation of the of the shell. We consider the case where the polytropic constant of the SN shell 
is y = 5/3 and we examine the evolution of a arbitrary shell perturbation. The dispersion relation is derived. The growth rate of 
the perturbation is found and its temporal and spatial evolution is discussed. The stability domain depends on the ejecta shell 
thickness, its acceleration, and the perturbation wavelength. 

1. Introduction 
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The Rayleigh-Taylor Instability (RTI) is a common phenomena in supernovae evolution. For a type II supernova (SN), the 
. neutrino-driven RTI is pr edicted during the cor e collapse. It is excited at the discontinuity surface betwee n the stalled shock w ave 
£NJ ' and the neutrino sphere djanka & Mtiller 1996I) . In supernovae interiors, for example for the SN 1987A ( Frvxel let al. 199 ll) . the 
■ RTI is thought to be responsible for the mixing of the heavy elements (Ni, Si) with the light ones (He, H). This mixing process 
( ^ ' allows one to obtain a better interpretation of the light curve of the SN 1987 A. 

. In this paper we consider the type II SN explosion. When the thermonuclear reactions cannot supply anymore the energy 
' required to maintain the pressure equilibrium in the star, its core collapses and produces a neutron star at its center. If the rotation 
• axis and the magnetic field lines are not parallel inside the neutron star, the former becomes a pulsar which gradually converts 
Q ' the spinning energy into a flux of electromagnetic waves and high energy particles. This pulsar wind blows up the inner layers 
" ^ . of the supernova and creates a dense shell expanding in the circumstella r medium - the SN remn ant (SNR). This structure of a 
shell expanding under the pressure of a central pulsar is called a plerion dWieler & Panagia"l9 80). 

This object experiences several hydrodynamic instabilities, especially, the inner shell surface is RT unstable during the ac- 
celeration phase. It is suggested that th e RTI at the inner ejecta surface is responsible for the shell fragmentation and for the 
filamentary structure in the Crab nebula dHester et al. 199fj|) . For type la and type I I SNe, the RTI arises a l so at the outer surface 
of an old SNR wh en the ejected shell is decelerated by the interstellar medium dVelazauez et al. 19981 IChevalier et al. 19921 
iHerantetal. 199 ll) . The plerion stability analysis is complicated for the pulsar pressure is non stationary, the density distribution 
in the shell is not uniform and, in addition, one has to take into account the velocity of the shell since it is in rap id expansion. 

Although the RTI has been studied for a long time since the pioneering work of Lord Rayleigh dRavleigh 19921) 
and Sir G. Taylor JTavi or 19501) . only a limited number of exact analytical solutions is known. How ever, the application 
of the RTI to t he SN explosion, has been considered earlier by Bernstein and Boo k dBernstein & Book 19781) . Reynolds 
and Chevalier dRevnolds & Chev alier 1984 and Blondin et al. dBlondin etal. 200ll) . Similar studies about the stability 
of target implosions h ave also been performed in the field of inertial confinement fusion dHattori et al. 198ol lHan 199 It 
iGoncharov et al. 2000l) (ICF). In these both domains the hydrodynamic instabilities develop in spherical geometry. Although 
the spatial and temporal scales as well as signs of the shell velocity and acceleration are different, the qualitative results of both 
series of publications are comparable. 

However, the stability analysis is not completed yet. The study by Hattori et al. dHattori et al. 19861) is limited to th e stagnation 
period of the shell evolution and it assumes a uniform density profile. Bernstein and Book dBern stein & Book 19781) carried out 
the study of the RTI for a shell in expansion, but they considered only an asymptotic limit of the RTI growth rate. Blondin et 
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al. dBlondin et al. 200 ll) considered only the instabilities that occur in the late phase of the SN evolution, when the ejecta shell 
is decelerated outside the pulsar nebula. The present paper provides another aspect of the study, i.e., the early stage of the SNR 
evolution when the shell is accelerated (and its thickness increases - see later) by the pulsar wind. On the other hand, this approach 
works also for the description in the compression of a collapsing shell. 

In this paper, Section|2provides the basic equations in the co-moving frame. Trough this procedure, we transform the non- 
stationary evolution of the shell into a steady motion. We perform the stability analysis for the case of an ideal mono-atomic gas. 
An analytical solution for the unperturbed spherically-symmetric flow is given in Section|3] We present the geometry and discuss 
about the physical parameters relevant to the conditions of a pulsar nebula-supernova remnant interaction. The stability analysis 
is performed in the co-moving coordinate frame and the dispersion relation is obtained in Section|4] Section[5]is devoted to the 
numerical study of the dispersion relation and the nature of the unstable solution is studied in details in section[6] The temporal 
evolution of the shell outer and inner surfaces are analyzed in Sectional The conclusions are given in Section[S] 



2. Zooming coordinates 

Let us consider a spatio-temporal evolution of the shell expanding under the central pulsar pressure. Our model equations consist 
in the continuity and Euler equations for the gas density p and the velocity v: 

d t p + V r • (pv) = 0, (1) 

d t v + (yV r )v = -~V r p, (2) 
P 

for an ideal gas with a polytrope equation of state p — K p y and where the subscripts "f" and "r" stand for the partial derivative. 
The differential operators in the spherical coordinates (r, 8, <p) read: 

1 , 1 
V r • (pv) = — d r (r pv r ) + — — d e (sin 6pv e ) 
r L r sin 8 

+ — — ~ fy(pv>). 
r sinfl 

V r p = {d r p, - d g p, — — d$p), 
r r sin 8 

Vg V(p 1 2 2 

(v ■ V r )v r = v r d r v r + — d e v r + — d^Vr (v g + V A 

r r sin 8 r v 

(V ■ V r )vg - V r d r Vg + — dgVg + ^— d^Vg 

r r sin 8 

1 , 

+-(v r ve-cot6lvi), 
r v 

(v • V r )v0 = VrdyV* + — dgVA, + — ^— cLv* 
r r sin 8 

+ -(V r V 4 , + COt8v g V 4> ), 

r 

where v r , Vg and are respectively the radial, azimuthal and tangential component of the velocity v. Even for a radially symmetric 
flow this system of equations is difficult to solve because of the non-stationary character of the evolution. 

E qs. dD and (13 ar e solved by using a time-depende nt transformat i on initially introd uced by Munier and 
Feix ([Munier &Feix 1983) and, thereafter by Bouquet et al. (Bouquet etal. 1 985: Blo ttiau et al. 1988) and nowadays referred 
in astrophysics as 'zooming coordinates' (Hanawa & Matsumo to 2000t iHennebelle 2001; Shadmehri & Ghanbari 200ll) . It is a 
specific coordinate transformation from the laboratory space to the co-moving, radially expanding frame. This new space is 
labelled with a hat symbol and it is defined by the following relations: 

r = C(0 f, 8 = 8, <f> = <f>, dt = A 2 (t) di, (3) 
p = D(t)p, p = B(t)p, (4) 

where A, B, C, D are the time dependent scaling functions; p and p are the density and pressure in the co-moving frame depend- 
ing on r, 8, (p and t. The velocity v in the new frame is given by the classical definition: 

dr 

dt 

According to Eq. the velocity transformation reads: 
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where the dot stands for the time derivative. These two reference frames coincide at t = t = 0, that is: A(0) = B(0) = C(0) = 
D{0) = 1. In the new frame the equation of state reads p = (D y /B)Kp Y and Eqs. Q and (0 become: 

d { p + V ? ■ (pv) = -0,p, (6) 
3?v + (v.V f )v = 4v^-n 2 v-V, (7) 

p T z 

where Qi, AT, Q2 and t are the source terms arose due to the non-inertial nature of the new reference frame. We have: 

J C D\ BA 4 
Qi = A 3 — + — ), A 7 = 



C D/' C 2 D' 
C A\ 1 CA 4 
C~A/' ~ ~C~' 

where Q,-(i = 1,2) has the dimension of a frequency, r is a time and A 7 is a dimensionless number. The scaling functions 
A, B, C, D can be found from invariance considerations and conservation laws. 

First of all, we are going to require that the EOS in the new frame remains invariant compared to the EOS in the initial 
frame. In other words, the fluid keep the same thermodynamical properties in the new space as in the physical one. The condition 
p = Kp y implies that D y /B = 1 and, therefore, B = D 7 . This is the first constraint for the scaling function. 

Second, requesting the same property for the continuity equation, i.e., the amount of matter is preserved from one space to the 
other, one deduces — and we obtain D = (1/C) 3 . The single equation left is the equation of motion. In order to keep invariant 
the pressure gradient force which originates from fundamental physics laws in the initial space, we take N = 1 . Combining the 
former three relations for the scaling functions, one has (R ibevre et al. 20051) : 

A = C (3y ~ 1)/4 , B = C~ 3y , D = C-\ (8) 

At this stage, we have expressed the three scaling functions A, B, and D in terms of the single one, C, which is the function 
governing the scaling between the radial coordinates r and r. However we have to satisfy two more conditions: the coefficients 
in the friction term and in the radial forces in Eq. Q become respectively: 

Oy-5)CC iiy - l)l2 =2Q 2 , (9) 
CC 3r ~ 2 = T~ 2 . (10) 

In order to have a stationary flow in the co-moving frame, both parameters r and Q2 should be constant. An inspection of Eqs. 
and J shows that the only way to satisfy these two constraints is to set Q2 = by choosing the polytropic constant y — 5/3 
which corresponds to a mono-atomic ideal gas. Then, the friction force (proportional to v) vanishes in Eq. Finally, the scaling 
function C(f) is obtained from the resolution of the differential equation (II Oi and A(f), B(t) and D(f) are derived from Eq. 
This is the simplest case, which allows an analytical study of the instability development for a non-stationary regime in the space 
(r, 0, 4>, t). 

For another value of y (y + 5/3), the scale function C is still defined by Eq. (llOi . but then the function ^2(0 is neither zero 
nor constant and it accounts for an effective friction force in the co-moving frame. However, this approach is not very fruitful 
since a time-dependant coefficient appears in the rescaled equations and the problem becomes more complex in the (r, f)-space 
than the initial one. 

This result is clearly in apposition with our procedure which consists in making more simple a initial complex problem. 



3. Unperturbed flow solution 



With this choice of y all coefficients in (|5Jl and are constants and one may look for a static, radially-symmetric solution for 
these equations in the co-moving frame. It can be found by setting v(r, f) = Vo = 0, p(f, t) = po(r), and p(f, t) = po(r). Then the 
solution reads: 



Po(r) = po(0) 
p (r) = Kpl(0) 



1 



?-2 \ 



1/ 



1/(7-1) 



and 



p2\7/(y-o 



(ID 



where the central density po(0) = [r 2 (y - l)/2/ryr 2 ] 1 ^ y [) . This solution is spatially bounded for r 2 > in a sphere of radius r\. 
To return to the laboratory frame one must find the scaling function C from Eq. ( 1101 . 
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In the case of a mono-ato mic gas where y _ 5/3, the solution to Eq. (1 1 Oi is fully analytic and it is of most interest for a wide 
class of astronomical objects (Bernstein & Boo k 1978llvan der Swaluw et al. 20041) . The scaling function C reads: 

c(t) = ^(i + ?) 2 + (12) 

where C(0) = 1 and the constant f3 = rC(0) characterizes the shell initial velocity vo(r) = /3r/r is the initial physical (r, t) space. 
The function C(t) is plotted in Fig.[T]for three initial velocities. In all three cases C approaches the ballistic motion [C(f) oc t] for 
t » r, but for t < r the behaviors are quite different. If the initial velocity is positive or zero, C(f) is monotonic and the flow is 
always in expansion. If the initial velocity is negative,/? < 0, there is a minimum radius, which corresponds to a stagnation at the 
time t r — t\/3\I(J3 2 + 1), followed by an expansion regime. 

Once the function C(t) is found, the integration of the fourth equation in provides a relation between the two times t and i : 

t = Tg(t) 

where g(t) = arctan^S + (J3 2 + 1) r/r] - arctan/?. (13) 

Figure|2]shows the relation between i and t for three values of /3. Although t varies from zero to infinity in the laboratory frame, 
t is bounded and it varies in the range [0, i max [ where the upper value t max depends on /3 and is given by i max = T(n/2 - arctan/J). 

As the function C is known, one obtains the other scaling functions A(f), B(t) and D(t) from Eqs. (|SJl. As a consequence, the 
non-stationary solution in the laboratory frame can be derived and we get: 



p(r,t)=p (0)C 



-3 



1 - 



p(r,t) = Kpl /3 (0)C- 5 



C(t) 2 



1 - 



2 > 


3/2 










r 2 




5/2 


cm]) 





(14) 



(15) 



v r (r,t) = -^- (/3 + (J3 2 + !)'-)■ (16) 

The solution describes a bubble expansion, with density and pressure decreasing with time while the radius increases. The 
velocity increases linearly within the bubble. The flow is compressible and its divergence, 

A- V 3C 3 (g + *M 

fey= c = T iWTM ' (17) 

decreases with time as 3/f when t — > oo, i.e., the vol ume of the flow grows with time. This solution could be found in a different 
way by imposing initially the radial velocity profile ( London & Rose n 19861) : v r oc r. However, the method we used here is more 
efficient to perform the perturbation analysis presented below. 

From the above solution one can construct a model for a supernova remnant blown-up by a pulsar wind. For that we remove 
the inner part of the bubble within the sphere of radius ?q < r\ in the co-moving frame. The pressure of the removed fluid 
is replaced by a radiative pressure, on the inner surface of the shell, attributed to the pulsar wind. The solution corresponding, 
therefore, to an expanding shell in the laboratory frame is shown in Fig. [3] The inner (resp. outer) radius is given by ro(f) = C(f) ro 
[resp. n(t) = C(f) hi 

The density, pressure and velocity profiles in the shell follow directly from Eqs. (114-1 - dl6l ). They are given in Fig. [3] The 
shell thickness L increases as L = C Lq where Lq — r\ - ro is the initial shell thickness. The pressure law at the inner interface of 
the shell describes the pressure produced by the pulsar wind: 



Pra {t) = Kpf (0)cr 5 



rim 512 



r 2 (0)) 



(18) 



Figure|3d) shows the time evolution of the pressure at the inner surface of the shell. In particular, for t » t the pressure decreases 
as p ro oc r 5 . This temporal dependence of the pressure can be related to the pulsar luminosity: L(t) oc -d(p r<j V)/dt, where V is 
the volume of the inner cavity bounded by the shell. Since the volume increases as C 3 , the derivative of the product p ro V varies 
according to C/C -3 . As a result, asymptotically for t » r, one finds that L(t) oc r 3 . This dependence is in good agreement 
with the classical model of pulsar luminosity. Indeed, according to Blondin et al. dBlondin et al. 200 lb. for t » r the pulsar 
luminosity decreases as where 2 < p < 3 and r is the characteristic life time of the pulsar JCamilo etaL2 000). Consequently, 
the pressure law which follows from our analytical solution ( II 81 is relevant for a class of pulsars experiencing a slowing down of 
their rotation velocity. 

In addition, the radial distribution of the flow velocity in our model is also in agreement with observ ations. For in- 
stance, the motion and structure of the filamentary envelope of the Crab nebula has been studied by Trimble (Trimbl e~1968t 
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ICadez et al. 20041) . The motion of filaments is largely radial and the velocity of each filament is approximately proportional to its 
distance from the expansion center, v °c r. On the other hand, these observations indicate that the filaments have been accelerated. 
It turns out that this behavior can be probably fitted with our analytical model by an appropriate adjustment of the parameters /3 
and t. 



4. Study of the shell stability 

In this section we study the linear stability of the expanding shell against Rayleigh-Taylor (RT) type perturbations, assuming 
that the shell is accelerated by the pulsar wind pressure. As a matter of fact, this configuration is RT unstable for the light fluid 
(radiation) which pushes the heavy fluid (shell of the supernova remnant). As it is shown in Fig. 0] small perturbations at the 
inner surface might grow during the acceleration phase and eventually break down the shell. 

The stability analysis of the shell in the laboratory frame is complicated since the background flow is non-stationary and the 
density and pressure vary with time and space. This analysis is much more simple in the co-moving frame since the unperturbed 
shell there is at rest. We define any perturbed physical quantity q(f, 9,<p,i) = qo(f) + Sq(f, 9, <p, t) where 6q{f, 9, <p, t) is the 
perturbation in the (r, f)-space, and the perturbed position of the inner boundary is 6f = fj(9, <p, t). Since Qi = and assuming 
y = 5/3, i.e., Q2 = 0, one then obtains the following set of linearized hydrodynamic equations [Eqs. (|6j-(0]: 



d f 6p + TjSKr po6v r ) + 



1 



r sin 6* 



<9g(sin 0p {) Sv e ) 



1 



; <9^(po<5v0) = 0, 



. - „ ,. (19) 

r sin 9 v 

po dj5v r = -d~ r 5p - T~ 2 r dp, (20) 

rp d{6vg = -dgSp, (21) 

f sin 9 po df5v,p - -d$5p, (22) 

where 6p = C 2 Q (f)5p and C 2 Q (f) = (5/3)Kp 2 ' 3 (r) = (f 2 - r 2 )/3r 2 is the square of the sound speed. 

One needs also to take into account the boundary conditions at the interfaces f^ an d ^1 ■ Let us assume that each of interfaces 
undergoes the perturbation 77, (i = 0, 1) . Then one should request the continuity of the pressure and the radial velocity at each 
interface, f = r ; + 77,, that is, 

d,-?7; = Sv r (fi), d-rpoih) f]i + Sp(n) = 0, 2 = 0, 1. (23) 

Introducing the relative density perturbation (density contrast) e = 6p/po, one can reduce the system of Eqs. (1191 1. (I20l >. fl2 li . 
and ( 1221 to a single second order partial differential equation (PDE): 



C 2 dh + 



2C 2 



If 
372 



d- r £ + 



2r 2 



3r 4 r 2 



e + 



C 2 



[sm6d g (smedge) + die]. 



(24) 



d;dv r = -C 2 d?e+ —e, 



r 2 sin 2 6 

In this PDE, if we keep only the first term of the LHS and the first term of the RHS, as well, one obtain the very simple PDE, 
<9?e = C 2 3 2 e, which just corresponds to the sound wave equation, as expected. In Ea.J24>. much more terms arise and they 
account for gravity dynamical effect existing in the initial (r, f)-space. 

From the equations of motion J20I >- (I22I> and from the relation between 6p and dp, the components of the flow velocity can be 
expressed in terms of the density contrast: 

If 

3T 2 ' 

^ 2 s0 d e e, (26) 
fsm6d f 6vt = -C 2 ^e. (27) 
The angular dependence of the perturbed quantities can be written in terms of spherical harmonics Yi,„(6, 0): 

e(f,0,cf>,t) = i(f,t)Y, m (0,<f>), (28) 
77,(0, <p, t) = Ht) Y lm (6, <p), i = 0, 1 , (29) 

where l(f, t) and £,(f); i = 0,1; are three functions to be determined later on, and where I and m are two integers with I > and 
m e [-1, 1]. Then Eq. J24l > reduces to the equation for the radial function %: 



fd;S\>g 



„,„ 2 s 



C 



.v() 



(2 



If 



3C 2 t 2 



8,-1 + 



If 2 



1(1+ 1) 



3t 4 C 4 



(30) 
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This is a linear PDE with time-independent coefficients. One is therefore allowed to consider an exponential time dependence 
for solutions. Denoting u the eigenvalue and S(f) the eigenfunction, we represent l(r, f) = (ao/fi)S(f)exp(a>i/T) and £,(f) = 
cii exp(<yf/r) (i - 0, 1) where a, is the amplitude of the initial perturbation at the inner (i=0) - outer (i=l) surface, and we assume 
that a, <k r,. Then Eq. ( I30i can be written in a dimensionless ordinary differential equation (ODE): 



d 2 S 12 _ 1R 
dR 2+ \R~ 1 -R 2 
(6(2R 2 - 1) 1(1 + 1) 



dS 
~dR 



\ (1 -R 2 ) 2 



R 2 



3co 2 
T^R 2 



o, 



(31) 



where R - fjf\ varies from Rq = ro/r\ to R\ = 1. This equation can be reduced to a canonical form with the transformation 
S(R~)=Ri*F(R 2 )/(l-R 2 y. 



t r^ 2p 
X(X - 1} d^ + 



x(fi + 3) - |- +fi 



dF 

dx 



U 2 + LL - 1(1 + 1) Li 2 +4ll- 1(1+1) + 3(0 2 \ 

a + a \F = 0, 

4x 4 1 



(32) 



where x - R 2 . It corresponds to the well-known hypergeometric equation provided the coefficient in the term F/x is zero, i.e, 
jj 2 + fi - 1(1 + 1) = 0. T here are two solutions for th e parameter yu: n\ - I and /j.2 = -(I + 1). Identifying Eq. (I32> with the 
hypergeometric equation JAbramovitz & Stegun 19721) . two linearly independent solutions write in terms of the hypergeometric 
function T(a,f3\ y, x) where the indices a, [5 and y are calculated for the two values fi^ found above. The general the solution to 
Eq. 13 H reads, therefore: 



S(R) = 



R< 



1-R 2 



d r(a u /3i;yi,R 2 ) 



■C 2 T(a 2 ,p 2 ;y2,R 2 ), 



1 -R 2 

where C\ and C 2 are two arbitrary constants, and where: 
ffi,2 = ~y + 1 + A, f3 h2 = '-j- + 1 - A, 



(33) 



n,2 



(34) 



with A = V/(/+l) + 4 - 3a> 2 . For the given set of indices, T reduces to the associated Legendre func- 
tions ( Abramovitz & Stegun 19721 . and we have T(a x ,fi x \y x ,R 2 ) cc R-'-^ 2 P(R)/ Vl - R 2 and T(a 2 ,{3 2 ;y 2 ,R 2 ) oc 
R MI2 Q(R)I Vl -R 2 where P(R) = P'^ 2 ( Vl - R 2 ) and Q(R) = Q 1 ^ 1 ^ Vl -R 2 ). As consequence, the solution 03 takes 
the simplified form: 



S(R) 



1 



[C, P(/?) + C 2 Q(R)] . 



(35) 



Vfl(l -/? 2 ) 3 ^ 2 

Next, the components of the perturbed velocity 6v(r, 9, <f>, t) = 6V(r, 9, (f>) exp (wf/r) can be expressed in terms of the function S : 



SV r = 

8V e = - 
8% = - 



3a>T 
a 



- ( l- R 2 ) ^L +2 RS 
dR 



3utR 
So 



3ojtR sin 9 



(\-R 2 )Sd e Y lm , 

(1 -R 2 )S8,Y lm . 



(36) 



Coming back to Eq. J35i to get the right S(R), the solution constants C\ and C 2 can be found from the boundary conditions [23] 
Dropping the exponential time dependence of the perturbations these boundary conditions provide: 



a S(R) = a 



3R 



R(l-R z ) 



1 -R 2 
dS_ 

~d~R 



R=R, 



S [2R 2 - (J(l - R 2 )], i = 0, 1. 



(37) 
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The derivative of S [see Eq. J35H can be found analytically ( Abram ovitz & Stegun 1973) : 
dS -S(R) 



dR R(l-R 2 ) 
A + l 



[(l-m(4-A)-3] 



(1 -R 2 ) 2 R^ 2 



[C 1 P 1 (R) + C 2 Qi(R)l 



where Pi (R) = P'^y 2 ( Vl - R 2 ) and Qi(R) = Q^ly 2 ( Vl - R 2 ). At the inner boundary, R = R , equations (03 provide two linear 
relations between C\ and C 2 : 



Ci P(R ) + C 2 Q(Ro) = 3Rl /2 yjl^R 



d Pi(R ) + C 2 Qi(Ro) 



3< 2 (1 



A + l 



to + A + 



2*5 



From these expressions one defines completely the constants Ci and C2 and one gets: 



3/?. 



3/2 



2(«o) 



A + / 



2 2 *o " 1 N 
a? + A + " 



-Gi(/fo) 



(38) 



3R 



3/2 



C 2 



P(«o) 



A + Z 



2 2«s-n 

w + A + - 



-Pi(Ro) 



(39) 



where S = Pi(Rq)Q(Ro) - P(Rq)Qi(Rq)- From Ref. dAbramovitz & Steeun 1973) . one can show that & is a constant equal to 
— r(A + /)/F(A - /) where T is the classical Euler gamma function. 

With C\ and C2 given by Eqs. fl38i and j39l >. the solution J35i is completely defined. It depends on the mode number /, on the 
position of the inner surface Rq and on the free parameter a>. The latter should be defined from the remaining boundary conditions 
at the outer interface R = R\. The perturbations of the motion both interfaces (R = Pq and R = Ri = 1) cannot be taken initially 
independently, and cannot further evolve independently. Indeed, from first equation fl37i at the outer surface R = Ri - 1 one 
finds: 



flj = -a lim[S(R)(l -R 2 )]. 

3 /?— *i 



(40) 



Notice that this relation does not forbid the calculation of the divergence of S near the outer boundary, because the density and 
the sound speed of the unperturbed flow go to zero there. Since in Eq. J35i . P and Q are rational functions of y — Vl - R 2 , one 
can expand S in Taylor series for y — > and Eq. J40l > becomes: 



3^ 
ao 



lim 

y-»0 



C 1 P(l) + C 2 Q(l) _ dP 
— + Ci — 



y 



d y y= 



y=0 



dy y=0 



(41) 



The first term in the right hand side of fl41i diverges if the numerator is non-zero. Therefore, the only way to satisfy this boundary 
condition and to obtain a finite displacement at the outer surface is setting the numerator to zero. This condition CiP(l)+C2Q(l) = 
corresponds to the dispersion relation for the shell instability we are looking for. Using the explicit formulae for C\ and C2, we 
obtain: 

2 A 2R l 

co = -A r 

l-R 2 . 



A + l P(l)Qi(Ro)-Q(l)Pi(Ro) 



[—2 P(R )Q(D-P(l)Q(Ro) 



(42) 
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It will be shown later in Sec.[6]that provided (1421 is satisfied, the second boundary condition at the outer edge is also immediately 
satisfied. 

Equation ( I42l i provides a relation between the growth rate of the instability, oj, and the mode number I. It is independent of 
the azimuthal number m, because the unperturbed system is spherically-symmetric and it will be shown later that, in addition, it 
does not depend on the shell inner radius Rq. 

Moreover, from Eq. J41i one obtains also a relation between the displacements at the shell surfaces: 

01 <+1 r[(A + o/2 + i] 

flo H(A-0/2] 
( [2 n(A + l) [jt n(A + l)) 

V n 2 \ 2 2 

This will be used later on to analyze the temporal evolution of the instability. 



5. Numerical solution of the dispersion relation 

The dispersion equation J42i looks rather complicated and we first consider its numerical solution. For this purpose, we represent 
it as a relation between two functions, Fi(A) = Fr(A), depending on the parameter A. Here, .FY (A) = 1(1 + 1) + 4 - A 2 = a> 2 is the 
left hand side of Eq. J4-2I > and Fr(A), its right hand side. In Fig.|5j Fi and Fr are plotted as a function of the variable A and we 
search for the intersection points of these two functions for given parameters I and Rq. In order to have an unstable solution (u 2 
should be real and positive) the intersection must occur in the upper right part of the plane (Fri, A) : < A < y/Kl + 1) + 4 = A,„ 
and < Fr < 1(1 + 1) + 4. This defines the upper limit for the growth rate for a given /: a> 2 < [1(1 + 1) + 4]/3. To have a stable 
solution, the intersection must occur in the lower right part of the plane (Fr^l, A) : < A < y/l(l + 1) + 4 = A m and Fr < 0. 
These restrictions are useful to constrain the numerical solution of the dispersion equation. 

An example of curves Fr, Fi is shown in Fig. [5] Both curves are surprisingly simple. The left function Fi is just a parabola 
which has a maximum at A = and becomes negative for A > A m . The right hand side has also a maximum at A = and it 
decreases with A. There are two intersection points: one corresponds to (J 2 > (unstable mode) and another one to w 2 < 
(oscillating mode at A = 6 for the present mode). The numerical solution indicates that both intersection points are independent 
of the shell thickness Ro. Especially, the point corresponding to unstable solutions is always located at A = 3 for chosen value 
I = 4. It will be shown below that this is a consequence of the nature of the perturbed flow. 



6. Behavior of the perturbed flow 

Until now, we have considered a general solution for the perturbed flow, without any assumption about its nature. Especially, we 
have been working with equations including compressible effects. However, this rather simple solution for the instability growth 
may suggest that special properties should probably arise. In this section we are going to study the properties of the perturbed 
flow. 

Since dg(sin 8dgYi m )/ sin 6 + (d^Yi,,,)/ sin 2 6 - -1(1 + l)T/ m , and after the injection of the second derivative of S, d 2 S/dR 2 
obtained from Eq. (13 11 . the divergence of the perturbed velocity 6V [given by Eqs. (I36H reads, 

,. ft flow f „ R IdS 2RS \) 

div<5V=^ — — > Y[ m . (44) 



hr \ co 2 \dR 1 -R 2 

It is worth to note that from the second boundary condition d37i the velocity divergence is zero at the boundaries. Now, we make 
a stronger assumption. Let us assume that the flow is incompressible everywhere, i.e., div 6V = for any value of the radius 
R. Then, from ( 1441 we can express the first derivative of S , dS /dR, in terms of S ans R and we can compute d 2 S /dR 2 , as well. 
Plugging these derivatives in the perturbation equation ( 13 11 1, we find after some simplifications that the coefficient in front of S 
does not depend on R and the differential equation reduces to a simple algebraic relation between u> and I : 

lo 4 -lu 2 -1(1+ 1) = 0. (45) 

This relation is a consistency condition between the eigenmode equation (13 1 i [it is automatically satisfied] and our additional 
constraint coming from (I44i . i.e., requiring that 6V = everywhere inside the shell. It is also consistent with the boundary 
conditions (1371 . Therefore Eq. ( 145 \ can be considered as the dispersion relation. We have to admit that up to now we have not 
been able to find a solution to Eq. (13 11 . satisfying (1371 . that differs from the one providing a zero divergence of the perturbed 
velocity. In other words the solution that is "constrained" by the boundary conditions (1371 describes an incompressible RT 
perturbation for which the dispersion relation is given by ( 1451 . 
The solutions of this quartic polynomial equation read easily : 



ft>l,2 



: V7TT, w 3 ,4 = ±i v7. (46) 
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Moreover, it should be emphasized that these solutions perfectly agree with the solutions of the general dispersion equation ( 1421 
found numerically in the previous section. As a matter of fact, for I = 4 we have found numerically a> 2 — 5 and co 2 = -4 leading 
to a> = ± V5 and a> = ±2i, respectively. These values are those exactly obtained from ( I46> by replacing I with the specific value 
1 = 4. 

As consequence, the solution we have found assuming a velocity profile with a divergence equal to zero does not look anymore 
as a "trick" to satisfy the whole set of constraints. This result is the proof that the RT unstable flow is indeed incompressible and 
the growth rate is unexpectively independent on the shell thickness. This result is quite new and is shown here in this paper 
for th e first time , although some people studied before in details thin and thick shells (see for instance the monography by 
Kull dKull 199 lb . All the modes are unstable for any value of the mode number I and the most unstable growth rate, among 
u>\(l) and u>2(l), is shown in Fig. [6] 

In addition, one can easily check that the flow is irrotational. This property arises since the unperturbed flow is radial. The 
vorticity of the perturbed flow in the laboratory frame l|5} is given by: curl 5v = curl 6v/C{t) since A(t) = C(t). As curl 5v = 
from Eqs. J36t . the perturbed flow is therefore also irrotational in the "physical" space. 



7. Evolution of the shell perturbations 

The effect of the instability can be further analyzed by observing the deformation f of the shell. It can be calculated from the veloc- 
ity field in the rescaled frame by using the relation dtf/di = Sv which provides £{f, 9,<f>,t) = {t/cl>)6v = (r/w) exp{u>i/T)6V{r, 9, (f>), 
where 6V is given by 03- 

Since the dispersion equation has four roots i46i . the general solution for the displacement is a linear 
superposition of four modes. In particular, for the radial displacement, one has : 



l(R, e, 4>, t) = a Yj e "'' /T Y '>» (0 > ®> 

where 



(47) 



!=1 



3<w 



-d 



dR 



are the displacement eigenmodes and A, are the mode amplitudes (Ai,2 are real, A3 = AT). Moreover, as the function S depends 
on oj 2 then £1 = £2 and £3 = £4. We have £,-(/?o) = 1 also. This equality comes from the normalization of with the inner 
initial displacement ao (see Section^). Then |/(1) is given by Eq. i43i . 

At this stage, we can calculate the deformation of the shell in the initial space (R, 8, <p, t). From the scaling relation for the 
distances [first equation of Eqs. ifj}] and from the relation J 1 3I > between t and f, the radial displacement g r (R, 9, <t>, writes : 



9, cf>, t) = a Q C(t)Y lm (9, <P) J] A/f«CR)e w ' gW , 



where ao = C(0)fio = &o- The radial velocity perturbation is provided by Eq.@ with A(t) = C(t) since y = 5/3 : 



5v = 



Sv 



6V 



C(f) C(t) 



(48) 



(49) 



where 6V is given by Eqs. d36l >. The derivation of d49l > requires further details. The velocities v(r, t) and v(r, f) in the (r, f) 
and (r, f) spaces, respectively, are connected through Eq. (|5j- In (r, f), we have v(r, t) = Sv, while the corresponding equation 
is v(r, f) = 6v + Vo(r, t) in (r, f) where Sv is the velocity perturbation and where Vo(r, f) is the background velocity given by 
Vo(r, t) = C(t)r/C(t). Introducing these expression in 0, one obtains d49l > immediately. In addition, we are going to introduce 
the quantity 6V(r, t) related to 5v by 5v = 6V exp[a>g(t)] (this relation is similar to the relation between 6v and 6V). 

It is now possible to derive the constants A, (/ = 1,4) from the initial velocity perturbation profile in the physical space. 
In the hat space, the radial velocity 5v r is given by 6v r = d£ r /di. On the other hand since we have 5V r (R, 0) = 6V,-(R, 0) at 
t = i = and for = = 0, the A,'s can be calculated from the boundary conditions taken as ^y(/?o, 0) = ao A), 0) = «o D\, 
SV r (Ro,0) = ao Vo, 6V,(1,Q) = oq V\ in the physical space where Dq,D\, Vq and V\ are four arbitrary parameters. Solving the 
corresponding system of equations, we find : 



Ai,2 
A 3 ,4 



1 



2(ii - b 2 ) 
1 

Wi - bi) 



Di - D b 2 : 



(Vi - V b 2 ) 



D bi -Bit i—(Yi - V b0 
VI 



(50) 
(51) 



where b\ = % 1,2(1) and b 2 = £3,4(1) and the eigenmodes are normalized by the condition |,(/?o) = 1 (i = 1,4). In the following we 
study the evolution of several specific perturbations. 
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First, we study the most unstable mode, oj\ = V/ + 1 by setting Ai = 1 and A2,3,4 = 0. The solid lines in Fig.0show the time 
evolution of the radial shell normalized displacement at the inner (a) and outer (b) interface, 770,1 /aoC(t), in the laboratory frame. 
It depends on two parameters: r - the characteristic expansion time and /3 - the initial shell velocity. The exponential law is valid 
only for t <K t (acceleration phase of the shell), indeed, exp[u>ig(t)] — > exp(wif/r), then 770,1/ao = ex P( f / T X with C(t) — > 1. For 
t » t, expg(f) — > exp[<x>i(7r/2 - arctan/3)] it following that 770,1 /ao(f/T) = exp[wi(7r/2 - arctan/3)] with C(t) — > f/r. The shell 
is in a ballistic motion and the perturbation increases linearly, i.e., 770,1 oc f/r - the function g(t) is bounded and C(t) is a linear 
function of time. In the figure[7]the displacement is divided by the scale function C(t) to exhibit clearly the perturbation growth 
due to the RTI. It is interesting to note that in the case of zero initial velocity the amplification factor e"" 1 ^ 2 does not depend on 
the expansion time. Moreover, the amplification is smaller if the shell is always in expansion. In opposition, the amplification can 
be stronger for initially converging shell (J3 < but | /3 |» 1). 

The numerical model/code employed here is called Pansy (M cCrorv et al. 19771) . It computes the time development of three- 
dimensional modes of coupled hydrodynamic, thermodynamic, and transport phenomena, including heat flow, viscosity, fully 
linearized about zeroth order spherically or cylindrically symmetric compressible flows. The zeroth order solutions are calcu- 
lated on a typical one-dimensional lagrangian grid and have the form / ; '(f), where the f s are all of the necessary hydrodynamic 
variables, and other variables including zone radius, and where j is the radial zone index. First order quantities, of the form 
f l J m (t)Y^ m {8)exp(im<p) for spherical geometry, are calculated with difference equations which are linearly perturbed forms of 
the former discretized zeroth order equations, rather than discretizations of the linearly perturbed continuous zeroth order equa- 
tions. This relatively conservative/Hamiltonian differencing approach produces considerably improved treatment of phenomena 
requiring high resolution, especially artifici al viscosity f or shocks in contrast with earlier form of Pansy which required higher 
resolution for the same accuracy ( Henderso n et al. 19741) . 

The analytical solution for the inner and outer shell displacement of the unstable mode has been compared with the simu- 
lations performed with the perturbation code Pansy. The agreement between the analytical solution and the simulations is very 
good. The difference is less than one percent. It is shown in Fig. that the inner interface growth is more rapid because this 
surface is RT unstable, while the outer interface is in contact with the vacuum and does not show a substantial growth. It looks 
stable. The deviation between the theory and simulations at the outer surface is larger than at the inner interface because it is 
more difficult for the code to handle the contact with vacuum. 

From Eq. (II 3i . the amplification factor exp[ti» g(f)] for t — > +00 is given by exp[ sjl + 1(tt/2 - arctan/3)]. For /3 = 0, we recover 
the special case studied by Bernstein and Book ( Bern stein & Book 19781) and the value they found is very close to the one we 
have. For example, for I = 40, Book and Bernstein calculated numerically the amplification ^ 10 4 , while our formula provides 
exp( VT+T7T/2) = 2 x 10 4 . 

The asymptotic amplification depends strongly on the initial velocity defined by the parameter /3 [see Eq. ( I16H . For example, 
if /3 — 1, the shell experiences a monotonic expansion phase, and the amplification decreases to 10 3 . Therefore, the shell is less 
fragile. However, for/3 = -1 the shell collapse initially and the amplification would be 10 s . As a result, the shell becomes much 
more fragile. We conclud e therefore that mo st dangerous stage of the shell evolution is the stagnation phase, which has been 
analyzed by Hattori et al. JHattori et al. 1986h . 

The purely growing eigenmode that we have studied above corresponds to correlated perturbations between both surfaces 
(inner and outer) which seems to be not too realistic. It would be more appropriate to consider either independent initial per- 
turbations of the inner and outer surfaces or the perturbation of the inner surface only. According to d50i and i5l\ . from four 
incompressible modes one can construct any kind of initial perturbation. Let us consider a case where the inner and outer inter- 
face displacements are opposite, i.e., Do = 1, D\ = — 1, Vo = Vi = 0. Hereafter we call this type of perturbation the "sausage". 
This "sausage" could be produced in the outer shell of the star by convection phenomena arising just before its explosion. 

From Eqs. and one finds: Ai = A 2 = (-1 - b 2 )/[2(bi - b 2 )] and A 3 = A 4 = (b x + l)/[2(&i - b 2 )]. The radial profile 
of the initial displacement within the shell is shown in Fig. (JSJl). The time evolution of such a shell perturbation in the laboratory 
frame is shown in Fig. for the mode / = 4. We have normalized the displacement by the scale function «o C(t) just to 
exhibit the amplification due to the RTI and suppress the shell thickness expansion. The growing perturbation at the inner surface 
of the shell due to the RTI is clearly shown, while the outer surface becomes spherical [see Fig. {8J))] because this interface is 
stable. The perturbation grows for time period satisfying t < r when the shell is accelerated and one can see that the displacement 
is larger along the polar axis 9 — [X axis on Fig. (HJ))]. 

Additional types of shell perturbations are interesting to study. Three of then exist : (i) the "kink" configuration for which the 
inner and outer initial perturbations are the same, Do = D\ = 1, Vo = Vi = 0; (ii) the "inner" shell perturbation for which only 
the inner interface is initially perturbed, Do = 1, Di = 0, Vo = Vi =0; (iii) the "velocity" alteration for which the positions of 
the inner and outer surfaces are initially unperturbed but just the inner velocity is perturbed, Do = Di = 0, Vo = 1, Vi = 0. This 
configuration describes, for instance, a wind pulsar flux fluctuation. 

The relative weight of the unstable mode in each of perturbation is given by the constant Ai. Hence, this parameter can 
be used to compare the configurations. For the "sausage" case with Rq = 0.5 and I = 4 shown in Fig. |8j Aj = 0.56. For the 
three other shell perturbations we find : Ai = 0.43 ("kink"), Ai = 0.50 ("inner"), and Ai = 0.22 for the "velocity" case. 
Therefore, all the three displacement perturbations configurations, "sausage", "kink" and "inner" are equally unstable. However, 
the amplification coefficient for these three modes is approximately twice smaller than for the growing single mode (Ai = 1 and 
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Ay = 0, j — 2, 4). These modes can probably be excited more naturally from the supernova progenitor heterogeneities. Finally, 
the "velocity" perturbation seems to be the less unstable. Its amplification coefficient is more than four time smaller. 

However, the respective importance of different perturbations depend on the shell thickness and the mode number. For exam- 
ple, for Rq = 0.9 and the same mode number / = 4 we find: A! = 1.41 for the "sausage" mode, A] = 0.22 ("kink"), Aj = 0.80 
("inner"), and Aj = 0.36 ("velocity"). In this case, the "sausage" mode is the most dangerous and it grows even stronger than the 
pure mode. 

By knowing the mode amplification coefficient, one can conclude about the shell fragility. In other words, we are able to derive 
a supernova shell fragmentation criterion due to the RTI. Assuming that the disruption of the shell occurs when the perturbation 
amplitude becomes comparable to the shell thickness, this criterion defines the critical amplitude, a cr , of the initial perturbation 
for the shell disruption. We have immediately this criterion in the hat space, with Eq. J47> for 6 — <p = 0, the displacement of the 
inner interface on the X axis for the more unstable mode becomes : $ r (Ro, 0, 0, f) = a cr Ai exp[wig(f)] = (h - ro) = h (1 - Ro). 
We can rewrite this criterion as : 

a cr 1 — Rq _ e (Av7+T ,;n 

— » e gv> v , t — > +oo. (52) 

h Ai 

This critical amplitude depends on the mode number /, the shell thickness, and the type of perturbation. The formula is valid for 
low modes, I < lu — 2nRo/(l - Ro), with a breaking occurring in the linear regime. The higher modes, / > Im, seem to be less 
dangerous. Although their amplification coefficient can be large, the corresponding perturbation enters in non-linear regime and 
grows less rapidly. In contrast, the lower modes, I < Im, are growing slowly and they lag behind the mode Im- Therefore, the 
linear theory predicts that the most dangerous perturbation corresponds to the mode Im and the critical initial amplitude is given 
by Eq. \52\ where / is replaced with Im- For example, for Rq = 0.65 and /3 — the most dangerous mode is Im ~ 12 with the 
"sausage" or "inner" type of perturbation (Aj w 0.5). Then we obtain a cr jr\ = 2.4 x 10~ 3 . The "velocity" mode is much less 
dangerous since it corresponds to Aj = 0.14. 

As result, a perturbation amplitude as small as 0.24% of the inner radius would be large enough to disrupt the shell. Figure|5] 
shows the initial shell with the critical perturbation a cr for a "inner" configuration (panel a) and the shell close to disruption at 
time t = 20t (panel b). The theory predicts that the shell will break down near the pole (X-axis). 



8. Discussion and conclusion 

In the first part of this paper, a non-s tationary spherical flow describing the expansion of SN ejecta has been derived by applying a 
rescaling method (Bouqu et et al. 198 5) to the Euler equations for a polytropic gas. This is a non-stationary solution that describes 
the motion of a shell with initial finite thickness Lq. The shell dynamics is described by two parameters: the initial velocity [3 and 
the expansion time t. By changing the sign of [3 one may study implosions and explosions. 

This solution is relevant for the description of the plerion evolution - type II SN remnant driven by a central pulsar pressure. 
The te mporal behavior of th e pressure in our model is in agreement with the classical spin-down power law of the pulsar lumi- 
nosity (Blon din etal. 200 lb . The parameter yS takes into account the initial kinetic energy of the shell, which is released during 
the explosion, and the time t describes the life time of the pulsar at the center of the ejecta. 

The transformation from the laboratory frame to the co-moving reference frame allowed us to perform an analytical study of 
the 3D linear stability of this time-dependant radial flow. The linearized hydrodynamic equations for a mono-atomic gas (y = 5 /3) 
have been solved using a decomposition of the perturbation in spherical harmonics. The dispersion relation defines the growth 
rate a> of any mode I in the co-moving frame, which corresponds to a finite amplification in the laboratory frame. The growth rate 
is independent of the shell thickness Ro and the azimuthal mode number m. Although no assumptions have been made concerning 
the perturbed flow, we found that the unstable perturbation is incompressible and irrotational. An analytical expressions for the 
growth rate and the mode structure were con firmed by comparison the analytical theory with numerical simulations performed 
with the perturbation code PANSY dMcCrorv et al. 19771) . The sign and the magnitude of the initial velocity play an important 
role in the RTI development. High initial velocities stabilize the shell. In opposition, for the case of initially collapsing shell, 
fi < 0, the perturbation grows to much higher amplitudes. 

By using a linear superposition of stable and unstable modes we studied the dependence of the amplification coefficient on the 
initial decomposition. It was found that the deformation for a thick shell of the inner interface is most dangerous and corresponds 
to a amplification coefficient twice smaller than for the unstable single mode. While for a thin shell the "sausage" configuration 
seems to be the more dangerous mode. For a given shell thickness Rq, we defined the mode Im which is the most dangerous for 
the shell disruption. A criterion is derived that defines the critical initial inner surface deformation that produces a fragmentation 
of the shell. 

A fil amentary structure in the Crab nebula expansion has been observed by Hester, dHester et al. 1 996) Sankrit et 
al. (iSankrit et al. 19981) . The structure of radial filament velocities indicates that they have a c ommon or igin - the SN spheri- 
cal shell. Such a filamentary structure has been reproduced from numerical simulations by Jun Jjun 199 8!). who shows that the 
RTI driven by the acceleration of the thin shell provides the main mechanism for the shell disruption. Our analytical model is in 
a good agreement with these observations and numerical simulations. By considering a shell with the aspect ratio 10 (Rq = 0.9), 
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the same as in the paper Jjun 1 998). we deduce that most dangerous mode is Im = 60. This agre es with the fastest growing mode 
observed in Jun's simulations. Moreover, the amplitude of initial perturbations (~ 1%) used in jjun 199 8) follows from Eq. \52\ 
for the realistic initial parameters: the velocity of shell is 500 km/s, the initial radius 0.2 pc, and the explosion time r = 500 yrs 
(this implies /3 = 1.3). This demonstration suggests that our analytical model is relevant to study the stability analysis of various 
SN remnants. 

Finally, this work can be applied to ICF or, more generally, to study of RTI in laser target design. 

In this respect, this study is directly us eful to Laborat ory Astrophysics issues and it can be used to design appropriate laser 
target to examine the RTI problem in SNR (Ribevre 2006). 
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Fig. 1. Scaling function C(t) in three cases: without initial velocity,/? = (1), with a positive initial velocity,/? = 0.5 (2) and with 
a negative initial velocity,/? = -0.5 (3). 




Fig. 2. Relation between the time i in the co-moving frame and the time t in the laboratory frame, for/? = (1), 0.5 (2), and -0.5 
(3). 
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Fig. 3. Density (a), pressure (b) and velocity (c) profiles versus the radius in the laboratory frame at t — (solid line), t = 0.5 t 
(dotted-line), and t — r (dashed-line), (d) plot of the pulsar pressure p ro at the inner surface of the shell versus time. The initial 
parameters are ft = and r ( j = 0.2 f\. The velocity is normalized by r\jj, the density - by the central density po(0) and the 
pressure - by Kpo(Q) 5 ^ 3 . 
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Fig. 4. Sketch of the shell explosion under the pulsar pressure. The inner, ro = ro C(f), and outer, r\ = f\ C(t), radii of the shell 
vary with time. The inner surface is pushed by the pulsar wind and it is RT unstable. 
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Fig. 5. Plot of the two functions Fl(A) and Fr(A, I, Rq) representing the left and the right hand sides of Eq. J42l >. The intersection 
of these curves provides the value of A, corresponding to the solution of the dispersion equation J42I . Both curves depend on the 
mode number / and Fr depends also on the shell thickness Ro (curves are drawn for Rq = 0. 1 , 0.5 and 0.8). However, we see that 
the intersection point does not move. In this example / = 4, the value of to 2 is 5 and -4, i.e., I + 1 and — / and we see that the 
intersection point is located at A = 3. 
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Fig. 6. Dependence of the dimensionless growth rate to on the mode number /, u = V/ + 1 [from Eq. d46ll. 
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Fig. 7. Comparison between the time evolution of the inner (a) and outer (b) displacements of the shell for the more unstable 
mode: A] = 1 and A2,3 ; 4 = in the case I = 4 with a inner radius Ro = 0.5 and (3 = 0. The displacement is divided by the scale 
function C(f) to exhibit clearly the perturbation growth due to the RTI. These plot show a good agreement between the analytical 
solution (solid line) and the simulation (dashed curve) with the perturbation code PANSY. The difference between the two results 
is less than one percent. 




Fig. 8. Plot of the initial "sausage" shell perturbation (a) at t — and (b) at t — r. Panel (c) shows the time evolution of the inner 
and outer shell displacements and (d) gives the initial displacement within the the shell. The mode number is I = 4 and the shell 
thickness is Rq = 0.5. 




Fig. 9. Plot of the initial inner shell perturbation at t = (a) and t = 20t (b). The mode number is / = 12 and the shell thickness 
is R {) = 0.65. The relative initial inner shell perturbation is 0.24%. 



